clear all
set more off 
set maxvar 15000 
clear matrix


	use "$Mydirectory1/3_Output/2_PooledData_analysis.dta", clear 
    keep if baseline_sample==1  

    estimates clear
    
    gen log_sq = log_father_baseline*log_father_baseline

    levelsof decade, local(decades)
    foreach d of local decades {
    eststo clear
    
        di "`d'"
        
    * Grab mean of income to evaluate
        sum log_father_baseline if decade==`d' [aw=wgt_sex_race], d
        local avg1 = `r(mean)'
        local r10 = `r(p10)'
        local r90 = `r(p90)'            
        
    * Regressions
        * Linear
        reg log_son_baseline log_father_baseline [aw=wgt_sex_race] if decade==`d', robust 
            local coeff = _b[log_father_baseline]
            local r_sq = e(r2)

        * Quadratic
        reg log_son_baseline log_father_baseline log_sq [aw=wgt_sex_race] if decade==`d', robust 
        estimates store log`d'
     
    * Other table components   
        foreach y in mean p10 p90 {
            if "`y'"=="mean" local cond "`avg1'"
            if "`y'"=="p10" local cond "`r10'"
            if "`y'"=="p90" local cond "`r90'"
            
            lincom _b[log_father_baseline]+(2*_b[log_sq]*`cond')
            estadd scalar eval_`y' = `r(estimate)' 
            
            local lb0 = `r(estimate)' - (1.96*`r(se)')
            local lb: display %-04.2fc `lb0'
            local ub0 = `r(estimate)' + (1.96*`r(se)')
            local ub: display %-04.2fc `ub0'
            
            estadd local bounds_`y' "(`lb', `ub')"  
        }
        
        estadd scalar r2q = e(r2)

        estadd scalar original = `coeff'
        estadd scalar r2o = `r_sq'
    
    }
        
    label var log_father_baseline "Linear coefficient"
    label var log_sq "Quadratic coefficient"       
        
    #delimit ; 
    esttab log1910 log1920 log1930 log1940 log1950 log1960 log1970 using "$Mydirectory2/appendix_d/table_main_estimates_log_quadratic.tex", replace     booktabs nostar label  keep(log_father_baseline log_sq) nobase noomit nonote alignment(c c) width(1.5\textwidth) se(3) b(3) brackets
    mtitles("1910s" "1920s" "1930s" "1940s" "1950s" "1960s" "1970s") 
    stats(original eval_mean bounds_mean eval_p10 bounds_p10 eval_p90 bounds_p90 r2o r2q N, 
    fmt(%04.3fc %04.3fc %-04.2fc %04.3fc %-04.2fc %04.3fc %-04.2fc %04.3fc %04.3fc %9.2gc) 
    label("Coeff., linear only"  "Slope at mean" " " "Slope, 10th perc." " " "Slope, 90th perc." " " "R-sq., linear only" "R-sq., quadratic" "Observations"));
    #delimit cr 
    
